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We present and numerically solve a new equation of motion for a charged particle under the 
influence of an external force, taking into account the radiation reaction. This covariant equation is 
integrodifferential, as Dirac-Rohrlich's, but has several technical improvements. First, the equation 
has the form of the second Newton law, with acceleration isolated on the left hand side, and the force 
depending only on positions and velocities: thus, the equation is linear in the highest derivative. 
Second, the total four-force is by construction perpendicular to the four-velocity. And third, if the 
external force vanishes for all future times, the total force and the acceleration automatically vanish 
at present time. We show the advantages of the new equation by solving it numerically for several 
examples of external force. 
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I. INTRODUCTION 

Accelerated particles carrying electric charge are a 
source of electromagnetic radiation. The momentum car- 
ried away by the radiation held affects the particles's clas- 
sical motion imparting a recoil force traditionally known 
as "radiation reaction" . It is somewhat frustrating that 
a century and a half after the first attempts, no fully sat- 
isfactory formulation of the motion of a charged particle 
exists that includes the effect of this back-reaction. Thus, 
theoretical efforts continue to this day pJ-Q. 

Most physicists are familiar with the Abrahams and 
Lorentz [6| equation of motion [23J 
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that incorporates a radiation reaction force ffs-x into 
Newton's second law, in addition to any external forces 
f ex t accelerating the particle. (Usually m and e will be 
the mass and charge of the electron, but the equation 
will be valid for any pointlike charged system.) Also well 
known is the relativistic generalization in terms of four 
vectors, known as Lorentz-Dirac equation 



mx^(s) = fW) + ^e 2 (x»{s) + (x(s)) 2 x»(s)) (2) 

whose derivation can be found in standard textbooks 0, 
H| . A few remarks will however be carried onto this work 
later on in section ITT1 

This equation has traditionally being recognized as af- 
fected by the problem of self-accelerated solutions, dis- 
cussed later in section Inland subsection IIVAI Therefore, 
it is often replaced by an integral equation, the Dirac- 



Rdhrlich equation 

x»(s) = J™ ds'eS s - s 'V L (^p- + (x(s')) 2 x^(s')\ • 

(3) 

This is an integro-differential equation for the accelera- 
tion, featuring a characteristic length, L = 2e 2 /(3mc 2 ) 
that for the electron is called the "electron relaxation 
length", L e = 1.876fm. The self accelerated solutions of 
Eq. © arc absent by construction (the formal integration 
is carried out with the condition that the acceleration 
does not increase exponentially at infinity). 

We present in this work a novel equation, equivalent 
to Eq.Q, but that has several important technical im- 
provements. This is 
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Equation ([3]) is indeed known to have several draw- 
backs. The first one is the apparent lack of causality 
as the acceleration depends on the force at future times, 
but this "advanced" formulation is not in itself a real 
difficulty as we will show in section Mil The second and 
most serious drawback is the phenomenon of "preaccel- 
eration" , at scales of order L where the particle acceler- 
ates before the external force begins. In this respect our 
equation will not improve the existing situation. Never- 
theless, our alternative equation (j4]) is superior in several 
important technical counts, that make its numerical im- 
plementation and solution much simpler. These are 

1. The new equation is linear in the highest (second) 
derivative, while Eq. ([3]) is quadratic in it. 

2. The acceleration x appears in both left and right 
hand side of Eq. ([3]), while in the new equation 
Eq. ((4]) the right hand side depends only on four- 
position and four- velocity, but the four-acceleration 
is formally solved for. 



2 



3. If the external force vanishes for all times, the new 
equation Q yields zero acceleration automatically. 
This happens also in Eq. ((3]) although it is not as 
obvious. 

4. The right hand side of Eq. ((3j) is orthogonal to the 
four velocity, but this is not evident. In the new 
equation ((4]), orthogonality follows trivially from 
the antisymmetry of the integrand under [i -<-> v . 
This orthogonality guarantees = and thus 
x 2 = constant at every step on an iterative solution 
without relying on cancellations not always obvi- 
ous. 

We will derive, analyze and numerically solve the new 
equation ((4]) for several simple cases of interest. Our 
numerical methods are describe in section Ivl and the ac- 
tual computations are reported in section PVTl Then sec- 
tion I VIII wraps up the discussion and summarizes our 
findings. 



II. DERIVATION OF THE NEW EQUATION 
A. Self-accelerated solutions 

The Lorentz-Dirac equation ^ features two radiation- 
related self-forces. One is mL(x(s)) 2 x fl (s), that ac- 
counts for the power radiated (see below eq. (|58p ) in 
an irreversible manner. The other, the Schott term 
mL'x fl (s), also present in the relativistic generalization 
of the Abrahams-Lorentz equation (see the derivation 
in [?| ) , represents a reversible positive or negative transfer 
of four-momentum between the charged particle and its 
near field. This structure entails that setting f^ xt {s) = 
for all times is not sufficient to guarantee x = in Eq. (|2|) . 
Moreover, the Schott term brings about self-accelerated 
solutions with exponentially diverging accelerations. 

Indeed, the "free" Dirac-Lorentz equation 



= L(x^( s ) + (x) 2 x») 



(5) 



can be multiplied by x^ , and using the orthogonality be- 
tween four- velocity and four-acceleration x • x = the 
radiated power term vanishes, leaving 
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that can be integrated to yield a family of solutions char- 
acterized by an arbitrary real constant b 



(7) 



This solution is self-accelerated, in the sense that for 
large proper time s the acceleration grows exponentially 
as e s / L , a clearly unphysical behavior. The integration 
constant b happens to be the component of the three- 
dimensional acceleration at s = that is parallel to the 
velocity, 



a(0) • v(0) = 6|v(0) | , as can be seen employing the re- 
duction in equation (|54p below. 

Thus, customarily the Dirac-Lorentz equation is com- 
plemented by a physical boundary condition requiring 
the acceleration to vanish at infinite time, lim s _ i . 00 |a| — > 
(free particle), or a less restrictive one asking that the 
product 



lim e- s/L \a\ 







(8) 



so that the aceleration does not grow exponentially. Any 
of these boundary conditions at infinity is sufficient to 
eliminate the self-accelerated solutions and yields the ex- 
pected behavior 



fUs)=0 Vs 



a5"(a) = Vs 



(9) 



B. Construction of the Dirac-Rohrlich equation 

Let us briefly remind the reader how the traditional 
equation of motion Eq. §3§ is obtained. Start by multi- 
plying Eq. ([2]) by an integrating factor e~ s l L . The Dirac- 
Lorentz equation then reads 



X 1 " \ = e 



_ -s/L ( fext 
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(10) 



The left hand side can be then written as a total deriva- 
tive 



ds V 



( fext 

\ mL 



(x) 2 x» 



(11) 



This can be integrated to yield an integral equation where 
the integration constant sets the acceleration at initial 
proper time Sj (since the starting point was the third 
order Lorentz-Dirac equation), yielding 



(12) 



ds ' e -s'/L ffex± + /-^ \ (s /) 
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In this expression we already see that the right-hand side 
builds up a four-force that is non-local in time. Imposing 
now the four conditions to the acceleration 
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i"(a) 







(13) 



that eliminate the self-accelerated solutions of the free 
Lorentz-Dirac equation, we can evaluate Eq. (|12j) at s = 
oo to read 

e~ s ^ L x^( Sl ) = J™ ds'e- s '' L (J^g- + (x) 2 x^ («') 

(14) 

that, upon substitution in Eq. (|12[) and combining the 
two integrals, results in Eq. ((3j). The self-accelerated so- 
lutions of the free equation satisfying Eq. (J7|) have been 
eliminated 24 1. This however is only certain when the 
external force is identically zero, for all times. 
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C. New treatment 

To obtain the new equation ((4]) we multiply Eq. (fTTj) 
by x"(s), which turns the vector family of equations with 
an index \x into a tensor family with /i, is, and we anti- 
symmctrize in the two indices fi, v. This manipulation 
has the merit of cancelling the non-linear radiation term 

Barut Q employed the same antisymmetrization in the 
absence of external forces to identify and then eliminate 
self-accelerating unphysical solutions, but did not extend 
it to the general case with external forces as we do here. 

We are left with a differential equation equivalent to 
the Lorentz-Dirac equation but where the non-linearity 
in the highest derivative has been traded for a higher- 
rank tensor structure 

±V{s) T s (' ' '''"'") - fa v ) = 

e s/L 

mL ^ ^ ext ~ ^ 1 ^ ext ^ ' 
The left hand side of this equation is an exact deriva- 
tive. Let us shorten the notation defining two auxiliary 
antisymmetric tensors, 

C"" = x fi x u -x v x fl (16) 
== f^V - f: xt x» (17) 
that, substituted in Eq. (|15p turn it into [25[ 



d_ 

ds 



{e~ s ' L C^{s) 



-s/L 
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-K^{s) 



(18) 



A formal integration leads to an expression analogous 
to Eq. (|12p above, in terms of the initial condition for 

We now restrict ourselves to forces that vanish when 
s — > oo and solutions satisfying the asymptotic condi- 
tion (that eliminates self-accelerated solutions of the free 
equation of motion), 



lim 



(e- s/L C^(s)\ 







so we obtain 

e- s/L C^{s) 
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ds'e~ s '/ L K^{s') 



(19) 



(20) 



To return this equation to one with the acceleration 
solved for, as in Newton's second law, we employ the 
property x iy C^ v = i M that follows from x 2 = 1. Multi- 
plying Eq. (|20j) by x v and applying this property, 



i"(a) 



mL 



(21) 



which is Eq. (|4]), and completes the demonstration. 

Finally we note the form that the equation takes if the 
external force is of electromagnetic nature, 

2" (a) = 



mL 



x v (s) / ds'e^' L (FZx p x» - FZx p ±n («') 



showing that the equation is cubic in the velocity (the 
position dependence will enter through the external elec- 
tromagnetic field). 



III. 



INTEGRAL FORMULATION OF THE 2 nd 
NEWTON LAW 



While the Lorentz-Dirac equation is local in proper 
time, eliminating the initial condition with a future 
boundary condition turns-on an integration over future 
proper time in Eq. §3§ and @. 

This however is by itself not an obstacle. Newton's sec- 
ond law itself can be cast in the very same form by means 
of one integration, so that the velocity is related to an 
integral over all future times of the force. In this subsec- 
tion we examine this formulation of the most basic and 
familiar law of mechanics, in order to dispell some causal- 
ity objections addressed here below and in section IIV Al 
and to test numerical iterative methods later on in sec- 
tion ESI 



A. Three-dimensional version 

Let us start this discussion by integrating the tradi- 
tional form of Newton's second law in three-dimensions 



m-v(t) = F(t) 
in an analogous form to equation (|12p. 



v(t) = v(ti) + - / dt'F(t') 



(23) 



(24) 



At this point the velocity is expressed as an integral over 
the force at prior times, back to the time chosen for the 
initial condition, so causality is manifest. 

We now select those solutions that are asymptotically 
free, so that the velocity takes a constant value at large 
time, limt^oo v(t) = v^. Evaluating equation (IM1) at 
t = oo we see that 



v(ii 



1 

m ., 



dt'F(t') 



(25) 



(22) 



that expresses the velocity at time U as an integral of 
the force over later times, in analogy to equations ([3]) 
and ([3]) (except these provide the acceleration and not 
the velocity) . This equation is obviously not in violation 
of causality, as one can reinterpret it easily as an initial 
value problem for Vqq. Thus, one could understand equa- 
tions such as ((3]) or (j4j as integral initial value problems 
for an acceptable (not exponentially increasing) acceler- 
ation at infinite future time. 

The whole point of this exercise is that advanced effects 
do disappear by rewriting Newton's second law in integral 
form in Eq. (|25[) as an equation for Vk, . If we had instead 
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imposed as asymptotic condition that the particle be left 
at rest for t = oo (v^ = 0), we would have obtained 

v(t) = -— f°° dt'F(t') (26) 
m Jt 

a somewhat mysterious advanced equation such as 
Eq. ([3]) or (TJJ in apparent violation of causality (which it 
is not). 

B. Four-dimensional version of Newton's equation 

Let us take now the relativistic generalization of New- 
ton's second law in Eq. (|2"3"|) 

mx^(s) = f^(s) . (27) 

This can be multiplied by an integrating factor ^e~ s / L 
(here L is any positive constant, not necessarily the elec- 
tron relaxation length) to yield 

- e - s/L f^{s) = e- s l L i^{s) (28) 

TO 

= 4 (e' s/L x t "(s)) + Ic-" /i i"(s) . 
as \ ) L 

We can group two terms into an effective force per unit 

mass, - ^ — x £* , and express the equation, integrating 

once with initial condition at Sj, as 

er s / L x"(s) = (29) 

again in total analogy with Eq. (|12[) . If we now re- 
strict the solution space to those whose Lorentz dilata- 
tion factor 7(1") = (y/T— v 2 /c^j cannot grow expo- 
nentially [26[, by imposing 

lim (e- s / L x"(s)) = , 

we can eliminate the initial condition by performing an 
evaluation of the integral equation (|29[) at infinity, and 
obtain 

^ w =-jr°°^-^(^-^) . (so 

Once more we write it as an advanced integral equation 
that has several similarities with our Eq. ©. First, in 
both equations the highest order derivative of the posi- 
tion [x in Eq. (f3"0)) . x in Eq. ([3])) appears on both sides 
of the equation. Second, the integrand in the right hand 
side of both equations has two terms, one depending on 
the external force and the other one on the movement of 
the particle (27j . 

But Eq. pOj) is different in one important count from 
Eq. ([3]): as we will show with several numerical exam- 
ples, the velocity solving Eq. (f30| remains constant at 
any time for which the external force vanishes, while the 
acceleration solving Eq. ^ begins changing before the 
force (preacceleration). 



IV. FURTHER THEORETICAL DETAILS 
A. Problems with causality 

The boundary condition in equation (TT31 and (TH)]) 
guarantees that the solution to the equation of motion 
does not grow up without bounds at infinity, and this 
eliminates the self-accelerated solutions from the free 
Lorentz-Dirac equation (and from the intcgrodiffercntial 
equations derived thereof), when the external force van- 
ishes for all times. 

Nevertheless, what happens if the external force f ext 
does not always vanish? 

The issue of preacceleration has been often remarked. 
Note from equations and (jH) that the acceleration 
at time s depends on the force at later times s' > s. 
Therefore, barring a built-in cancellation as happens in 
Newton's law in integral form, Eq. (j3T)]) . acceleration pre- 
dates the force: for a force vanishing before some given 
time, an observer that would have a time resolution of or- 
der L could in principle see the particle accelerate before 
its action. This will be visible in the numerical examples 
below. 

One often reads the argument that classical physics 
does not apply below the Compton scale, that is much 
bigger than the electron relaxation length 

he 

Acompton = — ^ 385fm > L ~ 1.876fm . 

TO e 

While true, two observations are in order. First, from 
a more practical point of view, the scale L in the equa- 
tion of motion divides proper time into units. From the 
point of view of the laboratory, the relevant scale is the 
dilated distance j(v)L. For a 3 GeV electron typical for 
example of synchrotron light sources (such as ALBA [ll[ 
in Barcelona, this factor is 5870, so that j(v)L = 0.011 
nm. This means that the preacceleration can be a phe- 
nomenon at the atomic, not nuclear scale. 

Second we should note that the electron relaxation 
length L has dimensions of length and a fixed value inde- 
pendently of H, while the Compton wavelength depends 
linearly on h. It would be puzzling that if h happened to 
be a factor 1000 smaller, preaccelerations would manifest 
themselves at length scales larger than the new Comp- 
ton wavelength. Quantum mechanical effects would not 
provide an escape in this case. 

Let us also briefly discuss the Dirac-Lorentz equation 
in the following form 

x" - Lx" L x 2 i^ . (31) 

TO 

If the radiation term in the right hand side was absent, 
this would be a linear inhomogeneous equation. Its gen- 
eral solution would have been an affine space spanned by 
the solutions of the homogeneous equation with f ext set 
to 0, added to any particular solution of the inhomoge- 
neous equation. The self-accelerated solutions described 
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in Eq. ([7]) are eliminated from the set of homogeneous 
solutions by setting the boundary condition at infinite 
time given by Eq.®. This is sufficient to obtain only 
physical solutions. 

However, when the non-linear radiation term is in- 
cluded, the boundary condition at infinite time is not 
sufficient, since, due to its non-linear nature, it couples 
the various modes of the linear equation. Thus, even if 
we have forced the self-accelerated modes to be absent 
near t = oo, there is no guarantee that they will not re- 
couple to the equation at finite time. They have only 
been eliminated from the free equation or from the linear 
equation without the radiation term. 

When solving the equation of motion numerically, one 
might come across an invariant interval s with f^ xt = 
but x M ^ 0. Then one knows that a self-accelerated mode 
is active. Since this is not a defect of our new equation, 
but of the entire theoretical set up, we will continue with 
this warning in mind. 

Some authors discussed in the past [lU, [l3| whether 
the preacceleration phenomenon is related to the lack 
of analyticity of the functions employed to expose the 
phenomenon (see for example fig. [H] below). The (short- 
ened) argument is that the Lorentz-Dirac equation re- 
quires certain analyticity hypothesis not satisfied by such 
forces (so the claim is that they should not be used at all 
with the Lorentz-Dirac equation). With analytic forces 
the particle at times before the beginning of the force 
is not really preaccelerated by a future force, instead it 
reacts to the analytic extension of the future-time force 
to present times. We are not totally satisfied by this ar- 
gument because Newton's equation in intcgrodiffcrential 
form, Eq. Q30p , has a similar mathematical structure and 
does not present this behavior, as seen in Fig. [2J 

We have constructed another clarificatory example be- 
low in Eq. (|56|) that presents the unwanted preaccelera- 
tion phenomenon while the force is exactly zero, as shown 
in figure |H1 The function is continuous and all of its 
derivatives are continuous on the real line, although it is 
not analytic in the complex plane. Still, it satisfies all 
conditions for the Lorentz-Dirac equation to be valid on 
the real line, thus weakening Valentini's argument about 
preacceleration being artificially induced by a poor math- 
ematical choice for the external force. 



B. Further integrations 

The class of non-linear integral equations of motion 
that we consider here requires an iterative solution for 
most example forces. Once we have performed one itera- 
tion of the integral in the relevant equation of motion, for 
example the new Eq. (j4|), we have the highest derivative 
in hand for all times, for example x(s). This function can 
in turn be immediately integrated to obtain the function 
and its lower derivatives, for which there is no incon- 
venient in setting initial conditions at initial time x(ti), 
x(U). 



Thus, the velocity at arbitrary time is calculated from 
the acceleration MS (X simple quadrature 

±"(«) = x"( Si ) + J ds'^(s') (32) 

and the position is obtained in an analogous manner once 
the velocity is known, 

x»(s) = x^(si) + J ds'x^is') . (33) 

Velocity and position will be updated iteratively at every 
step using these quadratures. 

All functions would be known at this point in terms 
of the invariant interval s. To plot them in terms of 
laboratory time, all one needs to do is keep track, for 
each s, of the pairs (x° (s) , x 11 (s)) , (x° (s) , x^ (s)) , and 
(x° (s) , x 1 * (s)) . The pair (s, x°(s)) itself provides the time- 
dependent Lorentz contraction factor 7 for all proper 
times and can be used to trade proper and laboratory 
time for each other. 



C. Limit of small L 

The limit L — > is physically most relevant as it corre- 
sponds to the non-relativistic limit (to — > 00 at fixed e). 
In this particular limit the equation accepts a perturba- 
tive solution complementary to the iterative one that we 
will pursue in subsection IV A[ 

Changing the variable under the integral sign from s' 
to a = (s' — s)/L, we have that the right hand side of the 
equation of motion m'x^ = F^(s) becomes 

/>oo 

F ft (s)=x v (s) dcre- J K fn '{s + La) (34) 
Jo 

with K^ v defined in Eq. (|T7|) . 
If we expand around L = 0, 

dK^(s) 

K^{s + La) ~ K^(s) + La . . . 

as 

that we substitute into Eq. (f3"4"|) , keeping only the zeroth 
order and first orders, and employing dae~ a = 1, this 
gives back 

mx^ ~ x v K» v + Lx v K» v (35) 

which is a form of the Lorentz-Dirac equation where the 
third derivative has been eliminated in terms of the ex- 
ternal force. Simple algebraic manipulations reduce it to 
the following form 

'»•'■" = fLt + I-''-. (fLti" - ./.",/ • ; ' / ' ) (36) 
- f^t + L{f^ t -(x-f ext )x^ (37) 

accurate to first order in L. 
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The first term returns the relativistic form of the sec- 
ond Newtonian law in the absence of radiation (e.g. neu- 
tral particle). Both terms are seen to be explicitly per- 
pendicular to the four-velocity, as befits a relativistic 
force dx^ ffj, = 0. (this property is kept order by order in 
the L-expansion) . Eq. (|36| is local; the acceleration at 
time s depends on the external force and its time deriva- 
tives and the velocity of the particle at the same time 
s. There is no preacceleration (it vanishes also in the 
intcgrodifferential formulation when L — > 0). 

If higher order terms in L are kept, one can express the 
radiation force on the particle in terms of the external 
force and its derivatives alone. 

In this same expansion, the radiated power can be writ- 
ten in terms of the external force alone as 



L 2 dp ext (s) 
m ds 



(38) 



(remember that the external force is a space-like four- 
vector, f 2 < 0). 



D. Limit of large L 

Also interesting is the limit of small mass m, so that 
in the Lorentz-Dirac equation ([2]) the inertial term mi 
is negligible against the radiation reaction terms (equiv- 
alently, L — > oo). In the Standard Model of particle 
physics (and thus, in all current empirical evidence) there 
is no charged particle with L larger than the electron's. 
However, this situation might arise in dimcnsionally- 
rcduced theories; for example, it is known that masslcss 
charged excitations arise in graphene sheets [H[). 

Now the Lorentz-Dirac equation 



becomes 



2e 2Jext 



{x 2 )x^ 



(39) 



(40) 



The tensors C^ v and K^ v can be defined as in Eq. (IT61) . 
and wc have 



d 3r 2 
_rv - _ k» v 

ds C ~ 2e 2K 



(41) 



Since no natural length scale appears, the exponential 
integrating factor damping the non-locality in time is 
now arbitrary, as in subsection MI Bl when treating New- 
ton's equation, so L is now unrelated to a relaxation 
length. Demanding weak asymptotic freedom at large 
s, e~ s l L C ilv — > 0, one can eliminate the future bound- 
ary term, and afterwards project over x v , as we did in 
subsection III CI obtaining the equivalent of Eq. (0| for a 
masslcss particle. The result is 

(42) 



Since the applicability of this equation to a physical sys- 
tem is not yet clear, we will leave its exploration for a 
future work. 



V. NUMERICAL METHOD 

A. Iterative approach 

The iteration of Newton's second law in integral 
form pop and also of the traditional equation of mo- 
tion ([3]) presents a difficulty not entangling the new equa- 
tion ©. 

To appreciate it, imagine applying the simplest ap- 
proach to the iterative solution of the equation of motion 
(we will eschew Minkowski indices and all other unneces- 
sary details, and refer simply to x as the solution to the 
dynamical problem, as in Newton's case, but this might 
be 'x if attempting a solution of the more technically in- 
volved Dirac-Rohrlich equation ([3])). 

In the naive-most approach, one simply guesses or cal- 
culates by other approximate means (e.g. the solution 
in the absence of forces, or the solution in the absence 
of radiation reaction, or an otherwise simple analytical 
case) a first guess x\^. Then one plugs this guess into 
the right hand side of the integral equation to obtain 



x m (s) = J ds'N(x [0] (s')) 

where N is the simplified form of 

N»(x(s')) = e^' L ( *W - 1^1) (43) 
\ L m J 

Iterating several times, the generic equation becomes 

x [k ](s) = f ds'Nix^is')) . (44) 



If, to within some criterion (ours is specified below 
in subsection IV Dp . convergence is reached, ||i[fc](s) — 
Z[k-i] i s ) 1 1 < e 7 one would have achieved a solution of the 
integral equation by iteration. 

However the appearance of the highest derivative on 
the right hand side of the older equation ([3]) brings about 
well-known instabilities, and this approach can fail. Con- 
sider the simpler force-free relativistic Newton's equa- 
tion, that should accept a solution with arbitrary con- 
stant velocity (set by the initial condition). The iteration 
of 



x(s) 



ds' 



e^ s - s '^ L x(s') 



(45) 



with a constant solution becomes 



^[fe-i] 



o(ss')/L 
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Analytically this equation is perfectly fine as the integral 
in the right hand side equals one. However, any finite- 
precision implementation on a computer means that the 
integral equals (1 — e) (with e either positive or negative 
depending on the numerical integration algorithm). 
Then one obtains 

X[k] = i[k-i] x (1 - e) 
that becomes, upon iteration, 

X[ k ] = x [0 ] x (1 - e) k . 

Depending on the sign of e, this converges (when the 
number of sweeps k is sufficiently large) to or diverges to 
oo, with all other constant values being unstable points of 
the system. Thus, valid solutions are lost upon iteration. 

Of course, this can be alleviated by Jacobi's method. 
To see it, let us recover the force and write the simpler 
rclativistic Newton's equation (|30|) in discretized form as 
a linear system 



Ax = b 



(46) 



with 



An, 



N 

£< 

n' = l 



,{s n s n ,)/L- 



-W„ 



L 



Here the integration weights have been denoted by w n i . 
Jacobi's method separates the diagonal piece with n = 



D n - 8 nn , L 

so that part of what would have been the right hand side 
in a plain iteration becomes part of the left hand side, 
yielding the vector equality 

Dx [k] =b-(A- D)x [k _ 1} (47) 

that is immediately solved by inverting D (a diagonal 
matrix). The convergence of this equation depends on 
the spectrum of D^ 1 . If the integration weights w n are 
small (tightly spaced grid) then S nn / dominates in D, so 
its eigenvalues are near 1 and no quick divergence desta- 
bilizes the code. Cancellations between the terms on the 
right hand side of Eq. (|4"7|) are then facilitated and fi- 
nite solutions can be found. Less favourable situations 
require more sophisticated algorithms. 

These difficulties are a consequence of the highest 
derivative appearing inside the integral in the equation 
of motion. The new equation (|4|) can be solved by simple 
iteration of the system and does not require the Jacobi 
method. However to stabilize the iteration we have found 
convenient to slow down the rithm of updates [l4[ , so that 
given xrfc_i], and after calculating the right hand side, 
J K(i[i c _i]), the update becomes, instead of Eq. ((33]), 

X[k) = a J A'(i[ fc _!]) + (1 - a)i[ fe _i] , (48) 
with the slow-down parameter a typically 0.5-0.75. 



B. Computation of integrals over invariant interval 

The grid over the invariant interval is composed of 
equidistant points s±, S2 ■ ■ ■ sjv- At the end of the compu- 
tation the grid is translated (for plotting purposes) into 
laboratory time t, and it ceases being equispaced given 
that the time-dilation factor j(v) is variable. 

To calculate the integral over the closed interval (s, oo), 
represented as (sj,...,Sn) on a discrete grid, we employ 
Simpson's rule (supplemented by a trapezoidal rule for 
the last integration interval if the number of points be- 
tween Sj and sjv is even). 

With finite computer time, s = oo must in practice be 
replaced by a cutoff s max on the maximum invariant in- 
terval of the particle. Since x(s ~ s max ) can no more be 
calculated from the integral equation (the later being ad- 
vanced), the choice of the cutoff has to be such that this 
calculation is no more necessary. This can be achieved, 
for example, if the external force vanishes after a certain 
value of s. Then one expects that the particle will see 
vanishing accelerations at later times, and the accelera- 
tion a can just be set to zero in the vicinity of the cutoff 
without need of calculation. Each physical case requires 
separate examination here to set the boundary condi- 
tions. This will be revisited in the next subsection IV CI 

The number of points of the grid has to be such that 
Si — s,_i < L ~ 1.876 fm, particularizing to an electron. 
This guarantees that the fast-changing exponential factor 
e 0-s )/ L i s integrated over with sufficient accuracy. 

If one wants to follow the motion of the electron at 
large, superatomic scales at the nanometer or above, one 
then needs several million points on the grid for a slow 
electron. However, radiation emission and radiation re- 
action are very small for slow electrons, and more inter- 
esting is the case of electrons with 7 factors of order 5-10 
(emitted by beta decay and other nuclear processes, or 
accelerated by small machines) up to the thousands (for 
synchrotron light sources) or hundreds of thousands (for 
colliders at the energy frontier). Since s max = t max /j 
for some average Lorentz contraction factor on the tra- 
jectory, to reach the nanometer scale in the laboratory 
the grid in s can be much coarser. We use steps of order 
0.1-1 fm in s and our calculations typically require 10 000 
points as an order of magnitude estimate. Such calcula- 
tions can run on a typical 2-3 GHz processor in few min- 
utes. Large computations running on a dedicated cluster 
or supercomputer can potentially reach human size scales 
even with fermi steps in invariant interval. 



C. Behavior at large times, initial condition 

Since the integral equation is built with a boundary 
condition at infinite time, it is worth delving on the be- 
havior of the computer code at large times. 

Vlasov [lH has observed that the solution with preac- 
celeration present and i M finite as s — ^ 00 is unstable to 
perturbations of the acceleration in the future. Our ex- 
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perience with the system confirms that setting and main- 
taining the boundary condition at infinity requires some 
care in the numerical analysis. 

Checking first Newton's Eq. ([30jl . it is easy to note 
that if we set f ext = at large times (asymptotically 
free-particle condition), then the zeroth component of 
the velocity satisfies 

±°{s)= f°° ds'^^-e^-^^ (49) 
J a L 

(for large s) and thus a constant solution is possible, and 
after pulling i° out of the integral, shifting this, and using 

f°° dr . T 

/ — e ' = constant = 1 , 
Jo L 

we obtain that i° = 7 (constant) at large time, so that 
the coordinate time scales linearly with invariant interval 
as expected. The same reasoning holds for the space-like 
components of Eq. (|30[) . yielding a constant velocity at 
infinite time, as expected from Newton's first law. 

To implement this behavior on the computer, all we do 
is require, after reaching a certain s f ree large enough, to 
read-off the last calculated point and prolong the solution 
continuously to the cutoff s ma x, so that 

X [S 6 (Sf ree , S max )) = i(sf ree ) . (50) 

In the next iteration, we start computing the integral 
right hand side of Eq. (|50|) for s < Sf ree , and at the end, 
set again the last points in the largc-s interval equal to 
the new value at Sf ree . 

The initial condition drifts from iteration to iteration 
in Jacobi's method and needs to be reimposed in the it- 
erative process. We do this by rescaling the entire func- 
tion, be it by employing a linear shift or a dilatation, 
depending on whether the initial condition was or not 
zero. The rescaling stops being necessary once conver- 
gence has been achieved, so that the integral equation 
and the initial condition have simultaneously been satis- 
fied. 

The advantages of our new equation (that accepts a 
simple iterative solution such as Eq. (|48|)) come again to 
the forefront; all we need to do to impose initial values 
of position and velocity (and zero acceleration in the ab- 
sence of forces), is to set a;ro](si) i[o](si), and i[ ](si) in 
the first iteration, and not recalculate them. The pro- 
gram keeps their value for its entire duration. 

D. Convergence criterion 

At each iteration [n] we need to decide whether the 
solution is converging (the integral equation at hand is 
closer to being satisfied) or not. As a diagnostic number 
we employ an integral measure for each of the position 
and velocity (for Newton's Eq. (|30p in relativistic form) 
and an additional one for the acceleration (for Eq. (jU). 



For example the convergence function for the velocity 
is, assuming Si = 0, 




(where we are not employing the Minkowski metric, but a 
sum over absolutely positive quantities) . In the limit n — > 
oo, one would expect, if the numeric method approaches 
convergence, that C± — > 0. In practice this is reached to 
very good accuracy after several dozens of iterations. 

We find that the iteration does not tend to drift be- 
tween solutions (a common problem in numerical analy- 
sis), but it may develop instabilities and loose all contact 
with an approximate solution of the integral equation. In 
such case, one must stop execution of the program in an 
orderly manner by testing whether C± is near 1 (since the 
criterion is constructed as a ratio of difference over sum, 
this means that the two functions if , , , and if , are not 

[n+l] _ [n] 

commensurate and one is dominating the difference in 
the numerator, which is a signal of numeric divergence). 
After the stop, the program is ran again with a different 
initial guess or a tighter grid. 



E. Parameter relaxation 

Since the equation of motion is non-linear in i, often 
a solution will not be easily found for a given force. This 
is for example the case for oscillatory forces that have 
complicated solution functions, such as the particle in a 
magnetic field that we study below. 

A workable strategy to find such solutions is to select 
a physical or a grid parameter and start with a trivial 
value. For example, extend the grid only as far as the 
first oscillation of the force, and for larger s cutoff the 
force. One can then increase the cutoff in small steps to 
include more and more of the force. 

Another example is the parameter L, in certain cases 
it is convenient to work with very large or very small L, 
and then relax it to its physical value. 

In either case, one lets the computer code find a fully 
converged solution for a given L or s max , then varies this 
parameter slightly, and find convergence again, proce- 
dure iterated as many times as needed until the neces- 
sary value of the parameter is reached. We have used this 
method extensively in the computations in section IVI1 



VI. NUMERICAL EXAMPLES 

In this section we present several numerical examples. 
We will use a system of units with h = c = 1 , and length 
measured in fm. Also F will represent the force per unit 
mass in the laboratory frame. 
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FIG. 1: Solution of Newton's equation in integro-differential 
form in three dimensions, Eq. (|26l) for a step force (dashed 
line). We plot the velocity (solid line). The numerical solu- 
tion coincides with the expected one imposing the velocity to 
vanish at infinity as boundary condition. (Arbitrary units of 
length). 



FIG. 2: Solution of the relativistic version of Newton's equa- 
tion in integro-differential form in four dimensions, Eq. (|30[1 
for a step force (dashed line) . We plot the velocity (solid line) . 
The numerical solution coincides again with the expected one, 
and in particular there is no preacceleration. (Arbitrary units 
of length) . 



Subsection IVI Al serves to gain experience with the al- 
gorithm with Newton's equation in integral form, be it 
the standard three-dimensional or the covariant form in 
Minkowski space. 

In subsection IVI Bl we will consider several examples 
with radiation reaction: in each one, we define a force 
F in the laboratory frame, iterate the covariant equa- 
tion (0} until the desired convergence is reached [28| and 
then we represent the results in the laboratory frame with 
the help of relations ([54]) and (f55j) . For simplicity we 
will consider first a step function, and then an oscillatory 
force. In all cases represented in this section, the force 
plotted is the force per unit mass. 

Subsection IVI CI demonstrates that the equation and 
algorithm can be scaled to larger systems much beyond 
the few-Fermi scale, showing their practical use for real 
problems. 

A. Solution to Newton's equation in 
integro-differential form 
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FIG. 3: Solution of the relativistic version of Newton's equa- 
tion in integro-differential form in Minkowski space, Eq. (|30[) 
for an oscillating force (dashed line). We plot the velocity 
(solid line). There is no preacceleration. In fact, Newton's 
first law is fully satisfied, the velocity being extremum at the 
zeroes of the force. (Arbitrary units of length). 



As a test of our numerical iteration algorithm, let us 
take as a first example the three-dimensional second law 
of Newton in equation ([26]) for a square well repulsive po- 
tential, depicted in figure[TJ The force (in the z-direction) 
can be written down in terms of two step functions as 

*; = -/(0(t-fo)-0(*-*i)) (52) 
and the general analytical solution is 

v z (t) = v - 

f [(t - t Q ) (0{t - t ) - 9(t - h)) + (h ~ t )9(t - h)] 

(53) 



(remember that we use the force per unit mass). 

This solution for the velocity (incidently, showing no 
preacceleration) for the particular values / = 1 and (as 
built-in in Eq. (|26|) ) imposing — 0, is computed itera- 
tively with Eq. (f2"fjf and plotted in the figure, coinciding 
with the exact analytical solution [2^] . 

Next we perform the same computation with the rela- 
tivistic Newton law, Eq. (|30j) . and plot the result in fig- 
ure [H Again causality is manifest and there is no preac- 
celeration. We set as initial condition v(U) = 0.9ce z and 
the force parallel to the velocity. 

Finally we solve again Eq. (|30l) with a periodic po- 
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tcntial truncated at two times to and t\, with the result 
visible in figure [3] The initial velocity is again parallel to 
e z and the oscillatory force is now perpendicular to this, 
along the OY axis. 

Again there is no preacceleration. Moreover, Newton's 
first law is manifest, as it is seen that every time for 
which the force vanishes features also an extremum of 
the velocity (vanishing acceleration). 

In conclusion of these first exercises, there is nothing 
unexpected in casting and solving Newton's law in either 
three-dimensional or Minkowski form as an integrodiffer- 
ential equation with an advanced boundary condition at 
t = oo. 



B. Elementary solutions with radiation reaction 

We now proceed to show example solutions of Eq. (@| ■ 
The components of the velocity and acceleration in the 
laboratory shown in the figures have been obtained from 
the four-dimensional covariant quantities that appear in 
that equation through the relations 



(54) 



(55) 



The calculation is reported in figure 3J The iteratively 
calculated solution of Eq. ((4]) shows the following. The 
particle begins accelerating before the force starts: in this 
preacceleration interval, energy is not conserved (there is 
no external work and nevertheless the particle increases 
its energy and at the same time emits radiation, although 
both contributions are very small); this violation cannot 
be experimentally measured due to Heisenberg's energy- 
time uncertainty relation. As the force starts, initially 
the acceleration is almost constant and approximately 
equal to the external for ce p er unit mass, the energy ra- 
diated being very small [301 ] . 

Further on, if the velocity in the Laboratory tends to 
one, the acceleration must necessarily go to zero; but if 
the external force stops before approaching this asymp- 
totic limit, as in figure @] the acceleration falls suddenly 
to zero before the force ends, what we could call "pre- 
breaking" , a phenomenon that we have not seen reported 
before (traditionally emphasis was put only on the preac- 
celeration). Finally, after the external force ends the ac- 
celeration vanishes, as it is evident looking to Eq. (£[]). 

In figure Owe show the preacceleration and prebreak- 
ing effects with more detail by looking closely to the times 
to and t\ where the force starts and ends its action. 

If classical electromagnetic theory could be applied to 
forces very intense, acting during very short times, the 
oddities would manifest themselves in a more extreme 
way: we show them in the following examples as an ex- 
treme manifestation of the mathematical properties of 
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FIG. 4: Solution to the new equation of motion with radia- 
tion reaction, Eq. @ for a feeble step force (top plot, dashed 
line: modulus of the three-force in the laboratory; solid line: 
acceleration). While the velocity is far from c (bottom plot), 
the acceleration is practically constant. It diminishes slowly 
as the particle becomes relativistic. We see distinctly preacel- 
eration. There is also a "prebreaking" before the force stops, 
but due to the logaritmic scale it is not visible in this fig- 
ure; we will show it in linear scale in fig. [S] (Units of length 
are fm; note that the OY axis represents inverse fm for the 
acceleration and for the force per unit mass.) 



the formalism, without any pretension of physical appli- 
cation. 

In figure [5] we employ again a step-force, affecting a 
particle initially at rest. The advanced features are ob- 
vious in that preacceleration is present before the force 
acts, but also the acceleration is dropping quite fast well 
before the action of the force stops. In the bottom panel, 
it even starts dropping before the force acts. This last 
effect is a consequence of relativity, as the particle's ve- 
locity in the laboratory approaches c = 1, the labora- 
tory acceleration needs to vanish. To show this we plot 
in figure [7] the relativistic kinetic energy per unit mass, 
T/m = j(v(t)) — l, where we see that the energy increases 
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FIG. 5: Detail of figure [4] where we show smaller time inter- 
vals around to and ti, to discern the effects of preacceleration 
and prebreaking. These effects are a consequence of the ad- 
vanced integration in Eq. Q. 



steadily. 

In figure[5]thc force is the same as in figure^ but we set 
as initial condition a particle propagating (parallel to the 
force) with initial velocity 0.9c, so the effective change of 
the velocity (shown in the bottom plot) is much smaller. 

We now address the objection of Valentini [l2| dis- 
cussed in subsection llV Al stating that the Lorentz-Dirac 
equation was initially deduced under the hypothesis that 
further derivatives of the acceleration existed, and thus 
one should not employ a discontinuous force. We show 
in figure [9] an example where the force is proportional to 
a function 



F(t) 



if 

e -(l fm 2 )/((t- 

if 



to)(*i-*)) 



t < t 

t e (t , 
t > ti 



(56) 



that, although not being analytic, is continuous and has 
continuous derivatives to any order throughout the entire 
t real line. In this extreme (very large force acting dur- 
ing a very short interval) the resulting acceleration and 




FIG. 6: Solution to the new equation of motion with radiation 
reaction, Eq. Q for a step force (dashed line: modulus of the 
three- force) . The acceleration (solid line) takes off at a time 
of order L = 2e 2 /(3m e c 2 ) before the actual application of the 
force, in violation of causality, as is well known from other 
formulations of the theory. The particle is initially at rest. 
The force in the bottom plot triples the intensity of the force 
in the upper plot. (Units as in figure 21) 



velocity (we now take vo = 0.001) are similar to figure [5] 
demonstrating that the discontinuity of the force plays 
no particular role: preacceleration and prebreaking are 
also present. 

In figure [TO] we solve the same equation 0$ with a pe- 
riodic force that is switched on and off (modulated) by a 
step function, so that only three periods of the sinusoidal 
force are active. The solution presents self-acceleration. 

The reader should remember that in the figures we plot 
the laboratory force F l = fl xt /j, that is not directly 
what enters the four-force in the equation of motion due 
to the extracted Lorentz 7(11) factor, so that the effec- 
tive acceleration is not sinusoidal, but depends on the 
velocity. 

The last of this series of exercises is reported in fig- 
ure [TT] where the electron relaxation length L is artifi- 
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FIG. 7: Kinetic energy per unit mass 7(v(t)) — 1 (solid line) 
for a square-hat force that acts on a finite interval of time, 
same as bottom panel of figure [6] 
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FIG. 8: Same as in figure[6l except the particle moves initially 
along the force axis with velocity v/c — 0.9. The actual accel- 
eration | a | (dotted line in the top plot, note the log scale) is 
thus much less pronounced. In addition we also plot the ve- 
locity (solid line in the bottom figure). (Units as in figure 
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FIG. 9: We soften the step-function force to make it contin- 
uous, and with all derivatives also continuous, by employing 
Eq. 1(55]). (Units as in figure [p. 
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FIG. 10: Solution to the equation of motion under a rapidly 
oscillating sinusoidal force of finite extent. The phenomenon 
of preacceleration and the asymptotic condition of vanishing 
acceleration at t — > oo are clearly visible. While the force is 
acting, the acceleration is dephased by about n/2 with respect 
to it, and they do not vanish simultaneously. (Units of length 
are fm). 



cially diminished by a factor of about 20 to approach 
the non-relativistic limit of subsection IIV CI The acceler- 
ation is now almost in phase with the force as expected 
(although at zero crossing of the force it has appreciable 
non- vanishing values). 



C. Solving the equation with radiation reaction for 
large systems 

In this subsection our aim is to demonstrate that the 
new equation (jj) including radiation reaction, in spite 
of having a natural scale of L = 1.876 fm for an electron, 
can be scaled to work with large systems. We will climb 
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FIG. 11: Same as figure [lOl but artificially reducing the ef- 
fect of radiation reaction by taking a much smaller L ~ 
O(0.1)fm instead of the actual value L — 1.876fm equiva- 
lent to 2a/(3m e ). As can be seen, now the acceleration is 
practically in phase with the force, and in agreement with 
Newton's equation for the motion in an electric field without 
radiation reaction, showing that the integral formulation and 
computer code are not responsible for the preacceleration in 
figure [lOl the cause is rather intrinsic to the theory starting 
with Eq. ([2]). Note also that the acceleration in that figure 
was damped with respect to this, as the radiation reaction 
much reduces the acceleration of the particle, the forces be- 
ing equal. (Units of length are fm). 



to six orders of magnitude larger lengths in laboratory 
time, to the nanometer scale. 

The external force in the next example will be caused 
here by an external electromagnetic field, thus 



•/ ex, 



(57) 



Our first example will again be an oscillating field, but 
its magnitude and physical size is loosely inspired by the 
relatively new concept of crystal undulator [171 ] . 

An undulator is a well known device employed in syn- 
chrotron light sources and free electron lasers to produce 
electromagnetic radiation. An alternating magnetic field 
with a period of a few centimeters forces electrons to 
wiggle and emit synchrotron light or X-rays, that are 
collected for applications. Our new equation with radia- 
tion reaction could be applied to these systems, but we 
will start here with a smaller device, for an exercise with 
less computing power, where in several thousand steps 
the effect of radiation reaction is clearly visible. 

In a crystal undulator, positrons (and somewhat less 
effectively, because of their attraction to the positively 
charged crystal sites, electrons) are channeled between 
crystal planes. The elementary charge can suffer a hard 
collision with another electron in the medium or with a 
nucleus. Such processes require a Quantum Electrody- 
namical description [TH, HH and are beyond our scope. 
However, the continuous interaction with the crystal os- 
cillating Coulomb field resulting in the emission of rela- 
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FIG. 12: Power radiated by an electron channeled along a toy 
crystal of physical nm dimension for energy 300 MeV (lower 
curve, red online) and 3 GeV (top curve, black online). The 
full Eq. © is used. 



tively low-frequency radiation (as opposed to hard pho- 
tons) is amenable to a classical description. 

We take as a model electric field E = Eq sin(27rx / a)e y 
with Eq ~ e/r 2 , and a typical nanometer length r ~ 
a = 1 10 6 fm. Thus the electron, with initial velocity 
voe x , perceives a perpendicular oscillating field. One can 
characterize this field by a length scale A# = j^h that is 
of order 2.8 x 10 10 fm. Thus, there is a hierarchy of scales 
LCacAfi and we are happy that the computer code 
can handle all three scales without special attention. 
Taking L as the scale of the computation, A^ is a huge 
number, so that the force oc is tiny, but if left to 
act for large times above a, it certainly has an important 
influence on the trajectory. 

The trajectory is not particularly interesting, but we 
follow the emission of radiation for electrons of two en- 
ergies, 300 MeV and 3 GeV respectively, in figure IT21 

To make the example more interesting, to this electric 
field we add an additional field of equal intensity but 
rotating in the following way, 

E' = Eq (sin(27rx/(10a))e a; + cos(27ra/(10a))e y ) , 

that is, with a period ten times larger. This doesn't have 
any claim of realism to describe an existing device, but 
again it is just used to demonstrate the capability of our 
equation. For real undulators the period should at least 
be a factor of 10 larger, approaching even the micrometer 
scale, but maybe our calculation is more directly appli- 
cable to so called "nanowigglers" [2(| that are conceived 
at the nanometer scale. 

At any rate, the motion of the electron emitting clas- 
sical radiation can be numerically followed, as we do in 
fig. [T31 concentrating again on the emission of radiation, 
that presents strong oscillations and interference effects 
between the two periodic fields. 

We now turn to our second example, inspired by the 
physics near a neutron star crust. We take a positron 
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FIG. 13: Power radiated by a 3 GeV electron channeled along 
a toy crystal with initial velocity along the OX axis. An os- 
cillating field along OY has period 1 nm. An additional field 
rotating in the XY plane, and of equal strength than the 
first, has period 10 nm (simulating a wavy crystal deforma- 
tion). Top: complete motion. Bottom: detail of early times 
(the periodicity of the motion can be appreciated). The total 
energy radiated is here much smaller than the total energy. 



(emitted in beta decay or by any other means) with an 
initial Lorentz factor of 7 = 10 to be injected in a con- 
stant and homogeneous magnetic field along the OZ axis 
of intensity Bo = 10 14 Gauss, at the upper limit of the 
conditions usually agreed for such stars. The character- 
istic length for the field intensity is Kb — -5- — 470 fm, 
so that the hierarchy of scales is not as marked as in the 
prior example. 

Neglecting first the radiation reaction, the trajectory 
is helicoidal around the magnetic field B = Be z , that 
projects to a circumference on the XY plane with radius 
r± = 4.7 x 10 3 fm. This trajectory provides a good zeroth 
order guess xiq\ for the program iteration. Of course, 
the energy and gyration radius will vary once radiation 
reaction is included. 

Figure [TJ] shows the velocity corresponding to hav- 
ing solved Eq. @ for this field. Initially we inject the 
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FIG. 14: A positron with energy 5.1 MeV and initial velocity 
along -^=(1,0,1) is injected into a region with a constant, 
homogeneous magnetic field B = Be z . We plot v x and v y 
(dashed lines, v y starting at zero), shifted by 7r/2, as well as 
the total v± (solid line) and v z (constant dotted line). 
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FIG. 15: For the same particle and field as in figure [T4l we 
plot the Lorentz contraction factor, or kinetic energy per unit 
mass (top plot). Since a magnetic field does not exert work 
on a charged particle, the energy lost must be carried away by 
the radiation. The computation of TZ is shown in the bottom 
plot. 
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positron with equal 




at the origin of coordinates (therefore the axis of the 
helicoidal motion is parallel to, but not coincident with 
the OZ axis). However, the emission of radiation slows 
down the electron so that v± decreases as seen in the 
figure (the maximum v x and v v values are decreasing). 

Although the total velocity remains above 0.9c (having 
started with 0.99c), we see in figure[15]that the energy per 
unit mass drops rapidly to less than half its initial value in 
only three oscillations, below the nanometer scale. This 
is natural since the radiated synchrotron power for a cir- 
cular trajectory (taken from Eq.(5.14) in Eq.O up to the 
metric's signature) from 

^=-f^ 2 (58) 
grows with the fourth power of the energy, 

so the trajectory will be much more severely distorted 
for large velocities. The loss of energy is compensated 
by this emitted radiation (bottom plot of figure [T5| since 
the magnetic field performs no work on the particle. 
The asymptotic condition in this calculation is imple- 
mented here simply by the fact that, as the particle ra- 
diates energy away, the Lorentz force x B vanishes 
with the velocity and the particle reduces its radius of 
gyration and is thereafter left to coast freely in the z di- 
rection with no external acceleration. So at some large 
Smax we simply set the acceleration to zero. 

Thus, we have shown in this subsection the possibilities 
that the new formulation opens to accommodate non- 
trivial calculations of physical interest. 



possibly a consequence of preacceleration, by which an 
electric charge can pre-increase its energy of motion by 
absorbing it from the near field, to shed it later, thus 
passing classically above nominally forbidden barriers. 
This concept has been demonstrated for the Lorentz- 
Dirac equation [22j and we arc curious about whether 
the integral formulation here presented also supports it. 
This interesting feature deserves a complete (and thus, 
labour intensive) analysis that we postpone for another 
work. Our treatment of this problem by means of the in- 
tegral formulation will avoid the conundrum of separat- 
ing self-accelerated solutions in the presence of a force, 
and will thus have different systematic uncertainties. 

We have found, confirming the work of other au- 
thors [22| . that discontinuities of the force function or 
its derivatives are a feature that is not essential for the 
counterintuitive behavior of the motion influenced by ra- 
diation reaction. 

The new equation (Q| is linear in the highest (second) 
derivative, that is explicitly solved for and not left im- 
plicit. Practitionners of numerical analysis in classical 
electrodynamics will no doubt find this and the other 
advantages of the new equation superior and useful. We 
look forward to the further applications of the theory. 



VII. SUMMARY 



In this work we have presented a new equation of mo- 
tion for a classical particle under the influence of radia- 
tion reaction, with technical advantages for its numerical 
implementation that may make it useful for systems as 
varied as the astrophysics of neutron stars, or the syn- 
chrotron light emitted by a crystal undulator or nanowig- 
gler. We have by far not exhausted the number of areas 
where the new equation could be of use. For example, 
intense laser fields are increasingly becoming available, 
and reaching 10 22 watt per square centimeter [2l|, they 
are believed to enter the regime where radiation reaction 
is important. Since the field is characterized by extreme 
oscillations, our worked examples could be extended to 
cover the motion of a charge inside such field. 

Another striking example where we should like to ap- 
ply our equation in the future is "classical tunnelling" , 
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